Chromosome-level genome and recombination map of the male buffalo

Abstract Background The swamp buffalo (Bubalus bubalis carabanesis) is an economically important livestock supplying milk, meat, leather, and draft power. Several female buffalo genomes have been available, but the lack of high-quality male genomes hinders studies on chromosome evolution, especially Y, as well as meiotic recombination. Results Here, a chromosome-level genome with a contig N50 of 72.2 Mb and a fine-scale recombination map of male buffalo were reported. We found that transposable elements (TEs) and structural variants (SVs) may contribute to buffalo evolution by influencing adjacent gene expression. We further found that the pseudoautosomal region (PAR) of the Y chromosome is subject to stronger purification selection. The meiotic recombination map showed that there were 2 obvious recombination hotspots on chromosome 8, and the genes around them were mainly related to tooth development, which may have helped to enhance the adaption of buffalo to inferior feed. Among several genomic features, TE density has the strongest correlation with recombination rates. Moreover, the TE subfamily, SINE/tRNA, is likely to play a role in driving recombination into SVs. Conclusions The male genome and sperm sequencing will facilitate the understanding of the buffalo genomic evolution and functional research.


Bac kgr ound
For sexuall y r epr oducing or ganisms, meiotic r ecombination plays a vital role in generating genetic diversity and ensuring segregation of homologous c hr omosomes. Recombination e v ents tend to be une v enl y distributed in man y species and fr equentl y occur in small genomic regions termed recombination hotspots [ 1 , 2 ]. Genomic c har acters like tr ansposable elements (TEs), GC contents, and PRDM9 binding are reported to be associated with recombination frequency and promote the formation of recombination hotspots [3][4][5]. Hotspots among mammals and e v en between relative species are poorly conserved, and crossover regions are fastevolving and possibly facilitate adaptive evolution [ 6 ]. T herefore , the study of recombination for each individual is necessary for the further functional and evolutionary research on animals.
The domestic water buffalo is an importantly economic animal resource . T he global population size of the buffalo is about 200 million, and they supply milk, meat, leather, and draft po w er in a gricultur al pr oduction for mor e than 2 billion people [ 7 , 8 ]. Water buffaloes feed the largest human population all over the world among domestic animals and are viewed as the most exploitative potential liv estoc k by the Food and Agricultur e Or ganization. Two kinds of water buffalo, including swamp buffalo ( Bubalus bubalis carabanesis ; NCBI:txid346063) and river buffalo ( Bubalus bubalis bubalis ), are classified. Swamp buffaloes are mainly distributed in China and Southeast Asian countries, serving as the primary draft animals for rice growing over thousands of years [ 9 ]. Their strong bodies ar e ca pable of enduring the heavy work in the field. Howe v er, high-quality food is often in short supply in its living environment [ 10 ], which may have contributed to the buffalo's higher digestibility of crude protein and fiber [ 11 , 12 ]. Along with the boost of a gricultur al mec hanization, buffaloes ar e optimized for meat or milk production [ 13 , 14 ]. Buffalo meat contains less fat and cholesterol in comparison with beef, suggesting that it can decrease the burden on the cardiovascular system and therefore increase the benefits to human health. Mor eov er, buffalo meat is effective for the treatment of diabetes described in the Chinese medical classic "The Compendium of Materia Medica" [ 15 ].
Although se v er al of female buffalo genomes have been finished [ 9 , 14 ], the genome of a male buffalo, including the Y c hr omosome, is absent. Genome assembly of the Y chromosome is a huge challenge because of its massive repeat content, half the sequencing depth due to its haploid nature, and high similarity with some re-gions of the X c hr omosome. Furthermor e, the absence of a male swamp buffalo genome hinders the detection of sperm meiotic recombination on the Y c hr omosome and the study of its influencing factors. To solve these problems, we sorted long reads from the Y c hr omosome by a computational method and assembled them separ atel y to gener ate a high-quality genome of the male swamp buffalo. We further sequenced 78 single sperms from the same male buffalo to provide the first whole-genome recombination map in buffalo. The high-quality genome, fine-scale recombination map, and subsequent analyses are likely to facilitate the genetic breeding of buffalo and promote the compar ativ e genomics r esearc h.

Genome assembl y, e v alua tion, and annota tion
Many mammalian genome projects prioritize sequencing female individuals (XX) over males (XY), as the ha ploid natur e of the Y c hr omosome r esults in half its sequencing depth. This can decrease the assembled contiguity and length of the Y c hr omosome [ 16 ]. Ad ditionally, the high n umber of re petiti ve sequences and the similarity to parts of the X c hr omosome make the Y genome assembl y mor e c hallenging. Recentl y, a computational method based on population datasets was de v eloped to sort long reads and generate genome sequences of the male-specific region of the Y c hr omosome (MSY) [ 17 ]. This method was applied to male buffalo, resulting in a total length of 9.3 Mb of buffalo MSY with an N50 value of 1.1 Mb. The remaining reads were further assembled, and all resulting contigs were polished with 170 × ( ∼450 G) short r eads. Compar ed to the pr e viousl y published buffalo genomes [ 9 , 14 ], our assembly exhibited the best continuity with a contig N50 of 72.2 Mb (Table 1 ).
We further sequenced ∼60 × Hi-C data to scaffold these contigs. Inter estingl y, a contig with a length of 7.6 Mb sho w ed a strong interaction signal with both X-and Y-contigs ( Supplementary Fig.  S1), which is assumed to be the pseudoautosomal region (PAR). The contig was phased by HapCUT2 using short-read, long-read, and Hi-C data. We aligned the 2 haplotypes onto the X c hr omosome of a female swamp buffalo [ 9 ] to determine their locations. Finall y, we gener ated a c hr omosome-le v el assembl y including 25 long pesudo-c hr omosomes (N50 = 120.0 Mb) ( Fig. 1 A,C and Supplementary Fig. S2). Among them, 8 c hr omosomes consist of onl y 1 contig (Fig. 1 A). Eight c hr omosomes contain telomeric repeats at one of their ends, and 2 autosomes (Chr3 and Chr5) contain telomeric repeats at both ends. We identified centromeric repeats in 16 c hr omosomes, and all of them are acrocentric except for c hr omosomes 1-5, whic h is consistent with karyotyping anal ysis [ 18 , 19 ]. Chromosomes 1-5 are homologous to 2 or 3 cattle chromosomes separ atel y [ 9 , 20 ], and centr omeric r epeats ar e located in all junctions. Based on the comparison with the female swamp buffalo genome [ 9 ], our genome closed 287 gaps (65.0 Mb, maximum length is 2.4 Mb) in the female genome (total 532 gaps) (Supplementary Fig. S3). Additionally, we found more transposons, especially LINEs (Long Interspersed Nuclear Elements), that reach se v er al kilobases in length and fewer unknown or other repeats in our assembly (Fig. 1 B). All of these results suggest the completeness of our genome assembly of male swamp buffalo.
We further estimated the completeness and accuracy of the final assembly and found that it captured 95.8% of the BUSCO orthologs (Table 1 ). Using Merqury [ 21 ], we obtained a quality value score of 41.3 for our genome assembly. We mapped the short reads of the transcriptome on the genome and found 98.3% of them could be aligned. The homozygous single-nucleotide polymor phism (SNP) r atio was a ppr oximatel y 3.39 × 10 −6 per base pair based on genomic short-read alignment. Besides, about 92% of the annotated Y genes in the bull genome (Btau_5.0.1) could be explicitly ( > 90% identity and > 95% cov er a ge) ma pped to the Y c hr omosome. To perform genome annotation, we combined 3 methods, including de novo , homology-based, and tr anscriptome-based pr ediction. In total, we predicted 22,608 protein-coding genes in the male buffalo genome (Supplementary Table S1).

Evolution of genomic elements
TEs are ubiquitous in eukaryotic genomes and play a fundamental r ole in sha ping genomic function and evolution [22][23][24][25][26][27]. In male swamp buffalo, TEs account for a ppr oximatel y half (49.39%) of the genome (Supplementary Table S2). Among them, the LINE/RTE-BovB subclass is the most abundant TE, with a proportion of 17.77%. LINE/RTE-BovB repeats in ruminants are believed to be tr ansferr ed horizontall y fr om r eptiles [ 28 , 29 ]. We inv estigated 6 ruminant species with high-quality genomes and found that swamp buffalo LINE/RTE-BovB repeats are more active recently in swamp buffalo than in other species (Fig. 2 A). The kim ur a v alue of LINE/RTE-BovB burst insertion is 0.03, and the corresponding time is about 1.36 Mya (Million years ago) under a mutation rate of 1.1 × 10 −8 per generation [ 30 ]. This burst time is close to the time when the 2 buffaloes (swamp and ri ver) di v er ged [ 9 ], indicating that it may promote the differentiation of the 2 buffalo species. We discov er ed that about 14,000 genes of swamp buffalo contained LINE/RTE-BovB repeats in their intronic regions, and LINE/RTE-Bo vB might be in volv ed in the r egulation of man y genes, whic h pr esumabl y contributed to the differentiation.
In addition to TEs, structural variants (SVs) offer an alternativ e a ppr oac h for genome e v olution b y influencing gene expression and phenotypes [31][32][33][34][35][36][37]. We mapped both swamp and river buffalo to the cattle r efer ence genome (ARS-UCD1.3) and used Assemblytics to detect SVs. We identified a similar number of SVs in both buffalo species (82,877 for swamp and 82,747 for river), of whic h 63,352 wer e shar ed and 19,525 and 19,395 were unique to swamp and ri ver buffalo, se parately. The total lengths of SVs were 160.74 Mb and 144.55 Mb in swamp and river buffalo, respectively. Apart from deletions , the a verage length of all other 5 SV categories (including insertions, repeat expansions, r epeat contr actions , tandem expansions , and tandem contractions) in swamp buffalo was longer than that in river buffalo (Fig. 2 B). To investigate the impact of SVs on genes in swamp buffalo, we studied the expression of genes with SV insertions across diverse tissues. We found that genes with SV insertions tended to have the highest expr ession le v els in the lung tissue ( P = 1.7E-05) (Fig. 2 C). We investigated the condition of swamp buffalo genes with unique SV insertions and still found the same trend ( Fig. 2 D). Our analysis indicates that SVs in swamp buffalo may have contributed to the de v elopment and evolution of the respiratory system.
The genome construction of the Y c hr omosome pr ovides an opportunity to study the evolution of the sex c hr omosome in buffalo. It has been reported that mammals' Y chromosome undergoes abundant gene conversion [ 38 ], which leads to sequence homogenization [ 39 ]. We illustrate the intrachromosomal similarities across the swamp Y c hr omosome in a circle map (Fig. 2 E). It is evident that the sex differentiation region (SDR) sequence is more homogeneous than that of the P AR. Furthermor e, we identified par alogous genes within the SDR and between PARs of the X and Y c hr omosomes and calculated the dN/dS value of these paralogs . T he dN/dS value in the PAR was  lo w er than that in the SDR (Fig. 2 F), indicating that the PAR was subjected to stronger purification selection against possible gene damage caused by homologous recombination between X and Y c hr omosomes.

Identification of recombination events and hotspots
To investigate the landscape of recombination events in buffalo, we sequenced 78 sperms from the same male buffalo with an av er a ge depth of ∼5 ×, in total ac hie ving 99.8% genome coverage. By employing a set of stringent filtering measurements and the donor's heterozygous SNP information, we identified a total of 1,934,008 high-confidence SNP loci. Using Hapi [ 40 ] softw are, w e inferred chromosome-level haplotypes and identify recombination spots for each sperm (Fig. 3 A). In total, we iden-tified 1,956 cr ossov ers with an av er a ge of 25.1 per sperm cell, which is similar to that in human studies [ 41 , 42 ]. Approximately 74.8%, 63.2%, and 42.1% of these cr ossov ers wer e arr anged into the interval of 200, 100, and 30 kb, r espectiv el y, indicating a high le v el of r esolution ( Supplementary Fig. S4). The distribution of distances between adjacent cr ossov ers was not uniform, with a peak at a ppr oximatel y 50 Mb (Supplementary Fig. S5). Compar ed to noncr ossov er r egions with a density of 66.3 PRDM9 binding motifs (CCnCCnTnnCCnC) per Mb, we found a higher density of 69.2 binding motifs per Mb ar ound cr ossov ers, indicating a potential role of PRDM9 in regulating meiotic recombination hotspots.
Recombination hotspots are crucial for ensuring the proper segregation of meiotic chromosomes and generating genetic diversity in offspring [43][44][45]. We calculated the recombination rate with a 3-Mb sliding window and identified 2 distinct recombination hotspots, both located on c hr omosome 8 (Fig. 3 B). These hotspot regions contained 31 genes. By performing functional enric hment anal yses in using DAVID [ 46 ], we found that the most significant functional term was biomineral tissue development ( P = 5.6E-4), which included 3 tooth-related genes (IBSP, SPP1, MEPE) (Supplementary Table S3). MEPE, in particular, is thought to be str ongl y positiv el y selected in herbivorous mammalian lineages and plays a crucial role in promoting the formation and mineralization of dentin, thus contributing to the strength of tooth structure [ 47 ]. Notably, buffaloes are known to efficiently utilize coarse feed, such as straw, sunflo w er cakes, and sprouts, and convert them into valuable animal products [ 10 , 48 , 49 ]. Recombination hotspots may provide genetic diversity to these tooth-related genes, but further experimental validation is required to confirm their functional roles.

Factors affecting the recombination rate
To determine which factor(s) have the greatest impact on recombination r ates, se v er al suc h as PRDM9 binding, TEs, and GC content have been investigated. We performed a correlation analysis between these genomic features and the recombination rates . T he effects of density and length wer e anal yzed separ atel y for genes and TEs. We found that gene density and length had almost equal corr elations with r ecombination r ates , but for TEs , the density was significantl y mor e corr elated than the length (Fig. 4 A,C and Supplementary Figs. S6 and S7). Ultimately, among the factors analyzed, TE density was identified as the most influential factor on r ecombination r ates in buffalo (Fig. 4 A-D).
Pr e vious studies have reported that TEs are also the main source of SVs [ 50 ]. Ther efor e, it is speculated that TEs may affect the formation of SVs by increasing the frequency of recombination. We further investigated the relationship between TE subfamilies and recombination rates as well as SVs and found that SINE/tRNA had a str ong corr elation with both recombination rates and SVs (Fig. 4 E,F). SINE/tRN A w as also found to be an important source of SV in pigs [ 51 ]. Ho w e v er, further e vidence is needed to validate the functional role of SINE/tRNA in both recombination and SVs of swamp buffalo.

Discussion
We present here the chromosome-scale genome of male buffalo, which exhibits better contiguity than published buffalo genomes [ 9 , 14 ]. In addition, we conducted whole-genome sequencing of 78 sperms from the same male buffalo and constructed the first r ecombination ma p for buffalo. The high-quality genome, particularly the Y chromosome, and the recombination map provide v aluable r esources for e volutionary, br eeding, and compar ativ e genomic r esearc hes of swamp buffalo. Our study could have significant implications for the a gricultur al sector, particularl y in r egions where swamp buffalo are an important liv estoc k r esource. Our r esearc h may also hav e br oader implications for the study of genome evolution and recombination, which can provide insights into the genetic mechanisms that drive species diversification and adaptation. The study has the potential to impact the dail y liv es of farmers thr ough its contributions to the br eeding of water buffaloes for meat and milk production. By identifying ge-netic v ariation r elated to desir able tr aits and using this information in breeding programs, farmers can improve the productivity and profitability of their herds.
The assembly of the Y c hr omosome pr esents a significant challenge due to abundant and lengthy r epeats, r educed sequencing depth, and high homology with some regions of the X chromosome [ 16 ]. In this study, we overcame these challenges by performing deep long-and ultr a-long-r ead sequencing ( ∼105 ×) for the male buffalo. We used the SRY software [ 17 ] to sort the long reads of the Y chromosome, and these reads were separately assembled to overcome the last factor. We identified the contig of the PAR through the interaction relationship of the Hi-C heatmap and phased them by combining the second-and third-generation reads and Hi-C data. Finally, we obtained the buffalo Y genome with a total length of 17.2 Mb, which is well mapped by 92% of the annotated genes in the bull Y genome (Btau_5.0.1). The assembl y pr ocess for the buffalo Y c hr omosome can also be a pplied to other animals and plants containing sex-specific c hr omosomes or fr a gments.
Meiotic recombination is well studied in model species [ 5 , 41 , 42 , 52 ] but less so in liv estoc k. We sequenced 78 buffalo sperms and identified 1,956 recombination events with an aver a ge of 25.1 cr ossov ers per sperm cell, which is similar to that of humans [ 52 ]. The fine-scale recombination map revealed 2 recombination hotspots on c hr omosome 8 with significantly higher r ecombination r ates than else wher e in the swamp buf- falo genome. Intriguingly, genes near these hotspots were most significantl y r elated to tooth quality. Given that buffalo's primary food source is low-quality food such as plant str aw, r ecombination hotspots may generate genetic diversities in toothassociated genes to better adapt to the consumption of crude-fiber diets.
Se v er al factors, suc h as PRDM9 binding, TEs, and GC content, can influence recombination rates. We found that TE density had the strongest correlation with the recombination rate of swamptype buffalo. Furthermore, SINE/tRNA, a TE subfamily, was found to have a significant effect on both recombination rate and SVs. We speculate that this SINE/tRNA subfamily may contribute to intraspecies or interspecies genetic variation by promoting recombination. Se v er al studies hav e shown that the ZnF domain of PRDM9 recognizes specific DNA motifs and is responsible for the formation of recombination hotspots [ 43 , 53-55 ]. Ho w ever, the rapid evolution of PRDM9 results in changes in the DNA sequence it binds to [ 56 ]. The 13-bp motif (CCnCCnTnnCCnC) in humans that we used may not be optimal for the swamp buffalo PRDM9 binding r equir ements, whic h could lead to a weaker effect of the PRDM9 binding sequences on recombination frequency than TEs. Further functional assays are needed to determine the binding motif of swamp buffalo PRDM9. Ne v ertheless, compar ed with other fac-tors except for PRDM9 binding, TE density has a r elativ el y high correlation with the recombination rate.
In the future, the genome and recombination map of male river buffalo could be constructed, providing insights into the divergent domestication features between the 2 subspecies of water buffalo and facilitating modern breeding for meat and milk production, as well as identifying genetic v ariation r elated to traits of interest. Additionally, further functional assays need to be performed to c har acterize the binding motif of swamp buffalo PRDM9, whic h may lead to a better understanding of the factors affecting recombination rates. We plan to continue investigating the genetic basis of important traits in swamp buffalo and to explore ways to use this information to impr ov e br eeding pr ogr ams and animal welfare. We also hope to de v elop ne w tec hnologies and methodologies for studying the genetics of nonmodel organisms.

Sample collection and sequencing
We sampled blood DNA from a local male buffalo in the Guangxi Zhuang Autonomous Region. To construct a high-quality genome of the male swamp buffalo, se v er al platforms, including Illumina, nanopor e, Bionano, and Hi-C, wer e used to gener ate a bulk of datasets. Bionano Saphyr technology was applied and DLE1 restriction enzyme was used for digestion. Illumina Hi-C technology was used in this study. For the construction of Hi-C libraries, the buffalo DN A w as digested with the restriction enzyme MboI and then was sequenced on an Illumina Novoseq 6000 platform ( RRID:SCR _ 016387 ) with PE100 reads. We generated about 466.1 Gb (174 ×) Illumina short reads, 271.9 Gb (102 ×) nanopore long reads, 561.4 Gb (210 ×) Bionano molecules, and 291.8 Gb (109 ×) Hi-C data (Supplementary Table S4). The Hi-C data were used to scaffold the primary genome assembly, and Bionano data were further used to manually check the order and orientation. The sperm was collected at the r epr oductiv e medical and genetic Center of the People's Hospital of Guangxi Zhuang Autonomous Region and sequenced according to the pr e vious study [ 52 ]. We also sampled 14 tissues , including dorsal muscle , lung, liver, spleen, tongue , kidney, heart, hindleg, foreleg, adipose tissue, conarium, hypothalamus, cerebellum, medulla oblongata, and  Supplementary Tables S5 and S6, separ atel y. The cortical divisions ar e in r efer ence to humans [ 58 ].

Separ a tion of long reads belonging to the Y chromosome
We selected short-read datasets from 59 male swamp buffaloes and 62 female swamp buffaloes from our previous buffalo population study [ 9 ]. The datasets and long reads of the r efer ence male buffalo were delivered to the SRY software (v1.5) [ 17 ] to identify Yspecific k -mers and separated long reads belonging to the Y chromosome.

Genome assembly
The long reads of the Y chromosome and other chromosomes of the male swamp buffalo were assembled with nextdenovo (v2.4.0) [ 59 ], r espectiv el y. All of the assembled contigs were polished by nextpolish (v1.3.1) [ 60 ] with settings ( −max_depth 270 for shortr ead ma pping options and −min_read_len 1k and −max_depth 200 for long-read mapping options) using short reads. We used juicer (v1.5.7) [ 61 ] to align Hi-C data onto the male buffalo genome and identified a PAR region candidate contig, ctg000160, that str ongl y inter acts with both X and Y sequences . T hen, the extrac-tHAIRS pr ogr am in Ha pCUT2 (-indels 1) [ 62 ] was used to phase the ctg000160 contig based on the alignments of genomic short r eads, nanopor e r eads, and Hi-C r eads. The 2 haplotypes were mapped to the X chromosome sequences of the female swamp buffalo using the mummer software [ 63 ], and the more similar one was considered to belong to the PAR of the X c hr omosome. Finally, we used 3d-dna (v180922) [ 64 ] with Hi-C data to anchor the contigs and manually adjust their orders in Juicebox as well as c hec k with Bionano data for generating a c hr omosomele v el genome . T he completeness and accuracy of the final assemblies were estimated using BUSCO ( RRID:SCR _ 015008 ) v5.4.3 [ 65 ], Merqury ( RRID:SCR _ 022964 ) v1.3 [ 21 ], and short-read alignment.

Detection of SVs
We utilized the nucmer pr ogr am in the Mummer pac ka ge ( RRID: SCR _ 018171 ) (v4.0.0beta2) [ 63 ] to perform genome alignments between male swamp buffalo (or river buffalo) and cattle . T he resulting delta file was deliv er ed to Assembl ytics (v1.2.1) [ 81 ] for calling SVs. We set the parameters of Assemblytics with "10000 50 1000000" corresponding to unique alignment length and minim um and maxim um size of SVs, r espectiv el y. We a pplied the chisq.test function in the R package for the gene expression comparison of SV-inserted or unique SV-inserted genes with all genes of the male s wamp buffalo. T he r esults ar e listed in Supplementary Tables S7 and S8. Notably, SVs in this study refer to fixed genomic differences betw een sw amp and river buffalo and cattle and not to variants within a population.

Calculating dN/dS
To compute the dN/dS value of genes on the Y chromosome, we used blastp (ncbi-blast-2.9.0 + ) with e-value < 1-E05 to generated protein alignments for genes in the PARs of the X and Y c hr omosomes as well as self-to-self alignments for genes in SDR. Optimal alignments other than to themselv es wer e consider ed as homologous gene pairs . T he yn00 in the PAML pac ka ge ( RRID:SCR _ 0 14932 ) (v4.9) [ 82 ] was further used to calculate dN/dS values of paralogs.

SNP calling
Sequencing short reads for each sperm were mapped onto the male buffalo genome using BWA (v0.7.17-r1188) [ 83 ]. Bam files for the same sample were merged using samtools ( RRID:SC R _ 002105 ) (v1.9) [ 84 ]. Duplicate r eads wer e r emov ed using the rmdup command in samtools with default parameters. We used samtools mpileup with settings (-C 50 -min-MQ 30 -min-BQ 30) to call SNPs of all 78 sperms together. Using the samtools mpileup and bcftools [ 84 ] filter command (-e "%QUAL < 30 || DP < 30 || DP > 200" -g 5 -G 5), we called SNPs for the male buffalo r efer ence . T he genotype of single sperm should be consistent with that of the paternal genome, so we selected heterozygous SNPs of the sperms consistent with the r efer ence heterozygous SNP site for the identification of recombination e v ents. To detect cr ossov er e v ents in PAR, we aligned both X and Y single sperm to the PAR of the Y c hr omosome to identify biallelic SNPs.

Identifying recombination events in sperm
To detect recombination events in sperm, the Hapi package [ 40 ] in R was used to process the heterozygous SNP results of sperm. We follo w ed the oper ations r ecommended by the Ha pi softw are step b y step. First, w e used the "ha piErr orFilter" function with default parameters to remove the potential genotyping errors of sperms. Second, heterozygous SNPs that were genotypes in at least 10 sperm ( n = 10) were selected for constructing the high-quality fr ame w ork b y the "hapiF rameSelection" function, separ atel y. Imputation of missing data was performed by the "ha piIm upte" function with settings (nSPT = 3, allo wN A = 0). Thir d, w e inferr ed and pr oofr ead dr aft ha plotypes by "ha piPhase" and "ha piCVCluster" functions. Multiple cr ossov ers (cv-links ≥2) within 1 Mb were filtered. We further adopted a maximum parsimony of recombination (MPR) strategy to eliminate incorrect crossovers by the "hapiBlockMPR" function. Fourth, c hr omosome-le v el ha plotype assembl y was ac hie v ed by the "hapiAssemble" function, and the haplotypes located at the end of the c hr omosome wer e polished using the "hapiAssem-bleEnd" function with default par ameters. Finall y, we identified cr ossov ers in sperm by the "ha piIdentifyCV" function based on haplotypes for each sperm. Notably, some recombination e v ents may not be accur atel y identified despite strict conditions for the process of sperm genotyping and recombination event identification.

Da ta Av ailability
The genomic sequencing r eads wer e deposited in the Genome Sequence Arc hiv e in the National Genomics Data Center, with the accession number CRA007045. Genomic (PRJNA907420) and transcriptomic (PRJNA907420) raw data are also available via the ENA. The genome assembly and gene annotation of the male swamp buffalo were deposited in Figshare [ 85 ]. All additional supporting data ar e av ailable in the GigaScience GigaDB database [ 86 ].

Additional Files
Supplementary Fig. S1. The interaction between the candidate contig ctg000160 in the PAR region and the contigs of the X and Y c hr omosomes. Supplementary Fig. S2. The heatmap (resolution: 500 kb) of the male buffalo genome . T he incr ease in inter action signal is r epr esented from y ello w to red. Supplementary Fig. S3. The distribution of closed gap lengths (bin size 5 kb). Supplementary Fig. S4. Distribution of the internal sizes. Supplementary Fig. S5. Distribution of distances between adjacent recombinations. Supplementary Fig. S6. Correlation between SV length and recombination rate. Supplementary Fig. S7. Correlation between gene length and recombination rate. Supplementary Table S1. Statistics of predicted protein-coding genes in the male buffalo genome. Supplementary Table S2. Analysis of transposable elements (TEs) in the male buffalo genome. Supplementary Table S3. Functional enrichment of genes around recombination hotspots.

Ethics Statement
Samples wer e pr ovided fr om collabor ators for r esearc h that was undertaken at Foshan University, permit FOSU2023001 from the School of Life Science and Engineering.